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Abstract: 

We present predictions for ttbb production at the LHC in next-to-leading order QCD. 
The precise description of this background process is a prerequisite to observe associated 
ttH production in the H — )■ bb decay channel and to directly measure the top-quark 
Yukawa coupling at the LHC. The leading-order cross section is extremely sensitive to scale 
variations. We observe that the traditional scale choice adopted in ATLAS simulations 
underestimates the ttbb background by a factor two and introduce a new dynamical scale 
that stabilizes the perturbative predictions. We study various kinematic distributions 
and observe that the corrections have little impact on their shapes if standard cuts are 
applied. In the regime of highly boosted Higgs bosons, which offers better perspectives to 
observe the ttH signal, we find significant distortions of the kinematic distributions. The 
one-loop amplitudes are computed using process-independent algebraic manipulations of 
Feynman diagrams and numerical tensor reduction. We find that this approach provides 
very high numerical stability and CPU efficiency. 
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1 Introduction 



The discovery of the Higgs boson and the measurement of its interactions with massive 
quarks and vector bosons represent a central goal of the ATLAS [HE] and CMS [3] exper- 
iments at the Large Hadron Collider (LHC). The present limits from direct searches and 
electroweak precision measurements favour a Higgs-mass range below the W-decay thresh- 
old. In this light-Higgs scenario, the Higgs predominantly decays into bottom quarks, and 
its observation is very challenging at the LHC. In the dominant Higgs-production chan- 
nel, i.e. in gluon-gluon fusion, the H — bb signal is completely obscured by a huge QCD 
background. 

Associated production mechanisms, where the Higgs boson is accompanied by a mas- 
sive gauge boson or top-quark pairs, feature more distinctive signatures that can be ex- 
ploited to reduce the background to the H — )• bb final state. These associated Higgs- 
production channels provide direct access to the interactions of the Higgs boson with 
gauge bosons and heavy quarks. Their observation would permit to test the electroweak 
symmetry-breaking mechanism. But the QCD background to associated WH, ZH, and 
ttH production followed by H — )■ bb decay remains a critical issue, which requires signifi- 
cant progress in two directions. The low signal-to-background ratios must be increased by 
means of improved selection strategies; and more precise descriptions of the background 
are needed in order to reduce systematic uncertainties. In this paper we address the latter 
issue by providing next-to-leading-order (NLO) QCD predictions for the irreducible QCD 
background to ttH(H — )■ bb) production. 

The strategies elaborated by ATLAS and CMS to identify the ttH(H — t- bb) signal 
[2H7] are based on the full reconstruction of the ttbb signature, starting from a final 
state with four b quarks and additional light jets. After imposing four b-taggings, a 
reconstruction of the top quarks is performed. This permits to identify two b quarks as 
top-decay products. The remaining two b quarks constitute a Higgs candidate, and their 
invariant-mass distribution is the relevant observable to find the Higgs signal. However, 
the presence of multiple b quarks and light jets in the final state represents a serious 
obstacle to the correct identification of the bb Higgs candidates. Realistic simulations 
indicate that only about 1/3 of the selected b-quark pairs have correct combinatorics, 
while the other Higgs candidates contain b jets from top decays or miss-tagged light 
jets. This so-called combinatorial background significantly dilutes the Higgs signal and 
increases its background contamination. The QCD processes pp — ?■ ttbb and ttjj are 
the main background components. The latest ATLAS and CMS simulations [2l|3], for 
30fb~^ and 60fb~^, respectively, anticipate a statistical significance around 2a (ignoring 
systematic uncertainties) and a fairly low signal-to-background ratio of order 1/10. This 
calls for better than 10% precision in the background description, a very demanding 
requirement both from the experimental and theoretical point of view. 

More recently, alternative search strategies based on the selection of highly boosted 
Higgs bosons, which decay into "fat jets" containing two b quarks, have opened new 
and very promising perspectives, both for \^H(H — )■ bb) [8] and ttH(H — )• bb) [9]. In 
the case of ttH, this novel approach might enable a better background suppression and 
increase the signal-to-background ratio beyond 1/3. Moreover, three b-taggings would be 
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sufficient to strongly suppress the ttjj contamination. In this case the background would 
be completely dominated by ttbb production. 

The calculation of the NLO QCD corrections to the process pp -> ttbb, ffist presented 
in Refs. [IOl[TT] and subsequently conffimed in Ref. [12], constitutes another important 
step towards the observability of the ttH(H — )■ bb) signal at the LHC The ATLAS sim- 
ulations of the ttbb background [2llll[5] are based on leading-order (LO) matrix elements 
and the scale choice /xr = /^f = -E-thr/S, where £^thr = 2mt + m^b is the partonic thresh- 
old energyj^ These predictions suffer from huge scale uncertainties: the LO ttbb cross 
section can vary up to a factor four if the renormalization and factorization scales are 
identified with different kinematic parameters [1]. However, in the case of the ttH sig- 
nal |iT3|IT¥j it was found that setting the QCD scale equal to half the threshold energy leads 
to fairly moderate NLO corrections {K ~ 1.2). The same behaviour was subsequently 
observed in two other processes that feature a final state similar to ttbb. At the scale 
/iR_F = -E'thr/2, the NLO QCD corrections to pp — )■ ttj [15] and ttZ [TB] at the LHC amount 
to ~ 1.0—1.15 (for pTjet ^ 20— 50GeV) and K ~ 1.35, respectively. On the basis of 
these observations one might have expected that LO ttbb predictions obtained with the 
same scale choice might have a decent level of precision, say at the 20-30% level. However, 
it turned out that the NLO corrections to ttbb production, for /i^ p = -^'thr/2, are much 
larger and correspond to a, K factor of about 1.8 [II1[I2]- Apart from the sizable impact 
on the ttH analysis, this big K factor suggests the presence of large logarithms that tend 
to spoil the convergence of the perturbative expansion. As we argue, this is mainly due 
to the fact that the scale /xr^f = -E'thr/2 does not provide an adequate description of the 
QCD dynamics of ttbb production. To cure this problem we introduce a new and more 
natural scale choice. This leads to a much smaller K factor and also reduces the residual 
scale dependence at NLO. Besides the issue of scale dependences, in this paper we also 
investigate NLO effects on various differential distributions and selection cuts that are 
relevant for the ttH analysis, both within the traditional ATLAS/CMS approach and in 
the boosted-Higgs framework. 

In addition to its phenomenological relevance, the calculation of the NLO corrections to 
pp — )■ ttbb constitutes also an important technical benchmark. The description of many- 
particle processes at NLO plays an important role for the LHC physics programme [TTpS] . 
Numerous Higgs and new-physics signals, as well as their background, are characterized 
by multi-particle final states. These processes often involve high powers of that give 
rise to very large scale uncertainties if NLO effects are not taken into account. The 
most prominent multi-particle reactions that require NLO predictions are summarised in 
the so-called Les Houches priority list [T71[18]. In the recent years, the technical chal- 
lenges raised by these calculations have triggered an impressive amount of conceptual 
and technical developments. In particular, in order to treat one-loop amplitudes with 
particle multiplicities higher than five, various methods based on tensor integrals [TUH2S] 
or on-shell reductions [SUHH] have been developed. Very recently, these different tech- 
niques have lead to the ffist NLO results for six-particle processes at the LHC, namely for 
pp — )■ ttbb [mill], the leading- [45] and the full-colour contributions [46] to pp — )■ Wjjj, 
and for the qq contribution to pp — )■ bbbb [17] . 

^The scale choice adopted in the CMS simulations [3l[6l[7] is not documented. 
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Figure 1: Sample tree diagrams contributing to the qq — )■ ttbb and gg — )■ ttbb channels. 




Figure 2: Sample pentagon and hexagon graphs contributing to qq — t- ttbb and gg — j- 
ttbb. The gg (qq) channel comprises in total 1003 (188) graphs, including 40 (8) hexagons 
and 114 (24) pentagons. 

To compute the virtual corrections to ttbb production we employ explicit diagram- 
matic representations of the one-loop amplitudes. A key feature of our approach is the 
factorization of colour structures at the level of individual diagrams. This permits to 
reduce the CPU cost of colour sums essentially to zero. Helicity-dependent structures are 
algebraically reduced to a common set of so-called standard matrix elements. In this way 
the sums over physical helicities are strongly boosted. Tensor loop integrals are related 
to scalar integrals by means of numerical algorithms that systematically avoid numeri- 
cal instabilities from inverse Gram determinants and other spurious singularities [201121]. 
The efficiency of the calculation is strongly increased by recycling a multitude of common 
subexpressions, which occur both inside individual diagrams and in tensor integrals of dif- 
ferent diagrams that share common sub-topologies. As demonstrated by the remarkably 
high CPU speed of the numerical code, these procedures strongly mitigate the factorial 
complexity that is inherent in Feynman diagrams. The real corrections are handled with 
the dipole subtraction method [501 - I53] , and the phase-space integration is performed with 
adaptive multi-channel methods [5H456] . Our results have been confirmed with the OPP 
method [ilM] and HELAC-ILOOP within the statistical Monte Carlo error of 

0.2% fT2]. 

The paper is organised as follows. Section |2] is devoted to technical aspects of the 
calculation of the virtual and real corrections. In Section [3] we present predictions for the 
LHC. In particular, we discuss the scale dependence and investigate NLO effects on the 
shape of several distributions. Our results are summarised in Section |H In App. |X] we 
outline the algebraic reduction of helicity structures, and in App. |B]we provide benchmark 
results for the matrix element squared in lowest order and including virtual corrections 
for one phase-space point. 

2 Details of the calculation 

In LO, the hadronic production of ttbb proceeds via the partonic processes qq — )■ ttbb 
and gg — )■ ttbb, which are described by 7 and 36 tree diagrams, respectively (see Figure [1]). 
The virtual NLO QCD corrections to these channels involve 188 and 1003 one-loop di- 
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Figure 3: Sample real-emission diagrams contributing to the channels qq — ttbbg, 
gg — )■ ttbbg, and gg — )• ttbbg. 

agrams, respectively. A few examples of pentagon and hexagon graphs are illustrated 
in Figure [21 The real emission contributions are induced by the partonic processes 
qq — )■ ttbbg, gg — j- ttbbg, gg — )■ ttbbg, and gg — )■ ttbbg. The gg channel involves 341 tree 
diagrams. The gg, gg and gg channels, which are related by crossing transformations, are 
described by 64 tree diagrams each (see Figure [3]). 

In the following we describe the calculation of the virtual and real NLO corrections. 
Each of these contributions has been worked out twice and independently, resulting in 
two completely independent computer codes. Top quarks are treated fully inclusively, 
i.e. we do not include top decays. Moreover we handle bottom quarks in the massless 
approximation, corresponding to the five-flavour scheme. However, we do not take into 
account the suppressed contribution from initial-state bottom quarks. 

2.1 Virtual corrections 

The virtual corrections are generated with two in-house Mathematica programs 
that reduce Feynman diagrams and generate Fortran?? code in a fully automatized 
way. One of the two programs relies on FormCalc 5.2 [57] for preliminary algebraic 
manipulations. Here we outline the underlying structure of the calculation, with em- 
phasis on colour /helicity structures and tensor integrals. In this respect, both programs 
are organised in a fairly similar way. Since the treatment of the gg channel is already 
documented in Ref. [10], we focus on the gg channel. 



Diagram-by-diagram approach 



The virtual corrections are obtained from the interference of the one-loop and LO 
matrix elements summed over external-state colours and helicities. This quantity is com- 
puted on a diagram-by-diagram basis. 



^^_y^{l-loop) (y^(LO)r 



col hel 



r L col hel 



(2.1) 



The contributions of individual loop diagrams (F) are evaluated by separate numerical 
routines and summed explicitly. The Feynman diagrams are generated with two indepen- 
dent tools, FeynArts 1.0 [58] and FeynArts 3.2 [59] . 
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Colour factorization 

One of the key features of the diagram-by-diagram approach is that the cost related 
to the large number of diagrams is compensated by the possibility to perform colour sums 
very efficiently. This is a consequence of colour factorization: individual (sub) diagrams 
consist of a single colour-stripped amplitude A^^^ multiplied by a trivial colour factor C^^\ 

^iv) ^ ^(r)c(r). (2.2) 

More precisely, each diagram gives rise to 3""^ colour-factorized contributions of type (12. 2p . 
where is the number of quartic gluon vertices in the diagram. These terms are handled 
as separate subdiagrams. However, most diagrams do not involve quartic couplings, and 
their colour structures factorize completely. For instance, the last diagram in Figure [2] 
involves a single colour structure 

^(hex) \ ^ jaibc ja2cd f^^^^\ (T'^Tj (2 3) 



b,c,d,e 



where ai,a2 and i3,i4,i^,iQ are the colour indices of the gg and ttbb external states, 
numbered in this order. All colour structures can be easily reduced to Kronecker symbols 
and Gell-Mann matrices T"" = A"/2 by using 

and other well-known SU(A^c) relations. In the gg channel, this reduction leads to a colour 
basis of 14 elements, 
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Ci4 = 6°''^°'^6i.^ig6irj^. 







(2.5) 

For Nc = 3, only 13 of these colour operators are independent owing to the relation 

6 12 

Ci4 = -2^C, + 2 5]C, + Ci3. (2.6) 

i=l i=7 

The summation over external colours is performed once and for all at the level of the 
colour basis and the LO matrix element. To this end, we compute the colour-interference 
matrix 

hi = X! ^kCi, (2.7) 

col 

and reducing the tree matrix element in colour space, 

>l(LO) = ^_y^(LO)^^^ ^2.8) 



5 



we build the interference of Ai^^'^^ with the elements of the colour basis as 

jcr = {M'-°'y = E4, {Mry . (2.9) 

col I 

Then, upon reduction of the factorized colour structure of the loop diagrams, 

M(n = ^(r)c(r) = ^(r) [^^ ^ip^^ ^ (2.io) 

we obtain the colour-summed loop-tree interference as 

EA^(^) {M^^o^y ^ ^(D (j2cPM^r'] , (2.11) 

col \ k / 

(T) 

where the coefficients are simple numbers. The colour-summed result is given by a 
combination of previously computed colour-Born interference terms fl2.9p . This requires a 
single evaluation of the non-trivial colour-stripped amplitude A^^^ of each (sub) diagram. 

Helicity structures are handled in a very similar way. The helicity-dependent parts 
of all diagrams are reduced to a common basis of so-called Standard Matrix Elements 
(SMEs), and helicity sums are performed once and for all at the level of the SMEs-Born 
interference (see below). The diagram- independent treatment of the helicity-dependent 
parts of loop graphs is made possible by the covariant decomposition of tensor integrals, 
i.e. by replacing loop momenta (in the numerator) with external momenta and metric 
tensors. 



Covariant decomposition and numerical reduction of tensor integrals 

Tensor one-loop integrals with propagators and P Lorentz indices are expressed in 
terms of totally symmetric covariant structures {g . . . gp . . . p}ji^^fp involving g^'^ and the 
external momenta Pi, ■ ■ ■ ,Pn-i, 

(2.12) 

with D denoting the number of space-time dimensions. For details of the notation we refer 
to Ref. [21]. To describe A^-point integrals with > 5, tensor structures with only four 
external momenta would be sufficient. However, in order to avoid potential instabilities 
due to inverse Gram determinants we use a redundant set of structures, including the 
metric tensor and N — 1 momenta. 

The virtual corrections to qq — )■ ttbb and gg — )■ ttbb involve tensor integrals up 
to rank P = 3 and P = 4, respectively. The one-loop amplitudes are expressed as 
linear combinations of tensor-integral coefficients Tj^ j^. The latter are evaluated by 
numerical libraries that recursively reduce them to master integrals using the methods of 
Refs. [2mi2lj^ . Avoiding an explicit reduction of analytic expressions to master integrals, 

^We note in passing that the reduction methods of Refs. [20l[2T] have also been used in the related 
calculation (62] of NLO QCD corrections to the 2 — >• 4 particle process 77 — >■ ttbb at a 77 collider. 
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this numerical approach prevents prohibitively large expressions and permits to adapt the 
reduction strategy to the specific numerical problems that appear in different phase-space 
regions. 

Tensor A^-point integrals with N > 5 are expressed in terms of lower-rank and lower- 
point integrals exploiting the four-dimensionality of space- time |20l|2T]| The tensor rank 
and the number of propagators are simultaneously reduced without introducing inverse 
Gram determinants. Consequently, the maximal power of inverse Gram determinants 
resulting from the entire reduction is given by the maximal rank of four-point integrals, 
which never exceeds four in renormalizable gauges. Scalar hexagons and pentagons are 
reduced to boxes using Melrose's method [60]. Tensor 4-point and 3-point integrals are 
reduced to scalar integrals with the Passarino-Veltman algorithm [61j as long as no small 
Gram determinant appears in the reduction. If small Gram determinants occur, alterna- 
tive schemes are applied [2T]E1 More precisely, we make use of expansions of the tensor 
coefficients about the limit of vanishing Gram determinants and possibly other kinemat- 
ical determinants. One- and two-point tensor integrals are obtained with numerically 
stable analytic expressions. 

Ultraviolet (UV) divergences are regularized dimensionally throughout, but infrared 
(IR) divergences are treated in different variants, which comprise pure dimensional reg- 
ularization with strictly massless light quarks and a hybrid scheme with small quark 
masses. The corresponding scalar integrals are evaluated using the methods and results 
of Refs. [63l|6l], and different regularization schemes are translated into each other as 
described in Ref. [65] . 

The calculation of tensor integrals is implemented in two independent Fortran li- 
braries. This permits to perform detailed cross checks, which confirm the excellent nu- 
merical stability of the reduction procedure. An automatic cache system is implemented 
that strongly boosts the reduction by recycling a multitude of tensor integrals among 
Feynman diagrams with common sub-topologies. The virtual corrections to gg — )■ ttbb 
comprise about 350 different scalar integrals, which require roughly 10 ms CPU time per 
phase-space point on a 3 GHz Intel Xeon processor. The calculation of the complete set 
of scalar and tensor integrals with and without cache system takes approximately 40 ms 
and 200 ms, respectively. 



Rational parts 

In D = A — 2e dimensions, UV-singular tensor integrals give rise to 1/ euv poles, 

tdN 

Tl,, = fl.,, + (2.13) 
Consequently, their D-dimensional coefficients need to be expanded in D — 4, 

fiml,. = mTi,, - (2.14) 



Rational contributions originate from D-dependent terms in tensor-reduction identities 



resulting in so-called rational terms that are proportional to the pole residues -Rj^ 



■^Similar reductions are described in Ref. [25] . 
Similar procedures based on numerical evaluations of specific one-loop integrals fl9 p25j or expansions 
in small determinants [24] have also been proposed by other authors. 
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and in the loop-momentum-independent part of the diagram numerators. The relevant 
expansions are automatically performed by means of a catalogue of residues of UV poles. 

Note that in fl2.14l) we have implicitly assumed that rational terms resulting from 1/e 
and poles of IR kind vanish. This is a non-trivial and general property of one-loop 
QCD amplitudes. More precisely, while rational terms of IR origin can be present in 
the wave-function renormalization factors, in truncated one-loop amplitudes they cancel. 
This holds within the 't Hooft-Feynman gauge and similar gauge fixings, as was explicitly 
proven in Appendix A of Ref . [TU] . 



Algebraic reduction of helicity structures and helicity sums 

The helicity structures encountered in the explicit evaluation of all Feynman diagrams 
are algebraically reduced to a common basis of SMEs as described below. The general 
form of SMEs for the gg — ?■ ttbb channel i^ 

Mrr. = Q^'''''-'''e^Me^,{p2) HPshpi ■ ■ ■ 7pMPi)] [HP5hp,+, ■ ■ ■ IpMPe)] , (2.15) 

where Q!^'^'^''^'"'^' consists of combinations of metric tensors and external momenta. These 
compact spinor chains permit to decouple helicity information from the remnant parts 
of the diagrams, so that helicity sums can be performed in a diagram-independent and 
efficient way. In practice, the colour-stripped part of each loop diagram [see fl2.10p ] is 
expressed as a linear combination of SMEs and tensor integrals, 

m 

N-1 

= E E f^^n-jpT^.,, + rational parts. (2.16) 

P ji,-,jp=0 

The coefficients )Cm\i...ip are rational functions of the kinematic invariants. These func- 
tions involve only denominators from intermediate-particle propagators and are free of 
spurious poles that might generate numerical instabilities. 

Helicity sums are performed at the level of the interference of the diagram-independent 
SMEs with the colour-projected Born amplitude (12. 9p . 

= E-^-M"^"^ = E^«E-^-(-^rT- (2.17) 

hel I hel 

This matrix is computed only once per phase-space point employing the Weyl-van der 
Waerden spinor formalism of Ref. [68] . Using Mkm one can directly obtain the colour- and 
helicity-summed contributions of each loop diagram in terms of its colour- and helicity- 
independent form factors J^^^ and the coefficients c^P of its factorized colour structure 

EE-^^'^ = T.^^' (EcPm,^ . (2.18) 

col hel m \ k / 



^For convenience we consider the crossed process ggttbb — 0, i.e. we treat all particles and momenta 
as incoming. 
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Owing to the high number of SMEs Xim and the complexity of the corresponding form 
factors J-'^\ the representation (12.181) yields fairly large expressions. For instance, the 
size of the numerical routines describing individual hexagon diagrams in the gg channel 
is of the order of 0.5—1 MB. The reduction of helicity structures to SMEs is one of the 
key aspects that determine the size and the speed of the code. In order to avoid possible 
numerical cancellations, this procedure is entirely based on algebraic identities that are 
free of denominators. The reduction algorithm consists of two main steps. The first step 
is based on process-independent identities in D dimensions. To reduce helicity structures 
of type (I2.15P we employ: momentum conservation; Dirac algebra; ordering of Dirac 
matrices inside Dirac chains; Dirac equation; transversality and gauge-fixing conditions 
for the gluon-polarization vectors, p^s^iipj) = for j = 1, 2. The basis of SMEs obtained 
with these identities contains more than thousand elements for the gg channel. 

After these manipulations in D dimensions we extract all rational terms performing 
the relevant expansions in D — 4. We then proceed with a second reduction step, based 
on four-dimensional relations. Specifically, we apply two alternative reduction algorithms 
that are based on relations derived from Chisholm's identity. 

The first algorithm is constructed along the lines of the reduction described in Ref. [10] 
for the qq channel. Each fermion chain is split into two contributions via insertion of chiral 
projectors Ci;± = (1 ± 7^)/2, 

[v{pz)^au{pA)\[v{p->)Ti,u{pQ)\= [v{p^)TaU}xu{pi)\[v{pr,)Ti,ujpu{p&)]. (2.19) 

\,p=± 

This permits to employ various relations of the type [T0|[66] 

7'^7"7^a;± ® 7^ = 7^a;± ® (7m7''7"w± + l''l^l^.^T) , (2-20) 

where the tensor product connects Dirac matrices that belong to different fermion chains. 
By means of such identities one can exchange Dirac matrices between chains that are 
connected by 7^^ ® 7^ contractions. As described in Ref. [10], using identities of type 
(I2.20p in combination with the above-mentioned D-dimensional relations (Dirac equation, 
etc.) one can obtain a rich variety of non-trivial reduction identities. In this way we 
have constructed a fairly sophisticated reduction algorithm (see App. that relates all 
helicity structures present in the gg channel to 502 SMEs. In spite of its efficiency, this 
reduction procedure has the disadvantage of depending on process-specific aspects, like 
the number of massive and massless fermion chains and the number of external momenta. 
It is thus important to investigate the trade-off between the obtained efficiency and the 
time-consuming task of designing the reduction on a process-by-process basis. To this end, 
we have implemented an alternative and much simpler reduction method. This procedure 
is entirely process-independent. It does not make use of chiral projectors, and consists of 
a single four-dimensional identity, 

^A'l^M2^At3^A'4^M5 _ ^MlM2^M3^M4^M5 _ ^AtlAt3^M2^M4^M5 _|_ g^J■l^J"^rJ^^2^^^3^^^5 _ ^MlM5^M2^M3^M4 
_|_ ^M2M3^Ml^M4^At5 _ ^M2M4^Atl^M3^M5 _|_ ^M2M5 ^A*! ^M3 ^M4 _|_ ^M3M4 ^Ml ^M2 

_ ^M3M5^M1^A'2^M4 _|_ gfJ.4f^5 _ ^gfJ-lt^2 gt^3t^4 _ ^MlM3^M2M4 _|_ gfJ-11^4 gfJ.2fJ.3^ ^MS 

_|_ ^^MlM2^M3M5 _ ^MlM3^M2M5 _|_ ^MlM5^M2M3'J ^M4 _ ^^MlM2^M4M5 _ ^MlM4^M2M5 
_|_(^MlM5^M2A(4^ ^M3 _j_ ^^MlM3^M4M5 _ ^MlM4^M3M5 _|_ ^MlMS ^M3M4 ^ ^M2 _ (^^M2M3^M4M5 
_g^^.2^^^ g^^3^^5 _j_ ^M2M5^M3M4^ ^Ml ^2 21) 
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which can be derived from Chisholm's identity and permits to ehminate any spinor chain 
involving more than three Dirac matrices without introducing 75 and e-tensors|§ In this 
case we could reduce all gg-channel helicity structures to 970 SMEs. 

Comparing the number of SMEs obtained with the pro cess- dependent and process- 
independent algorithms, we observe that the former is superior by roughly a factor two. 
Thus, if we naively assume that the CPU efficiency scales with the number of SMEs, 
we would expect a factor- two difference in the speed of the numerical code. In contrast, 
we find that the CPU efficiency obtained with the two reductions is almost identical. 
This suggests that the reduction of the number of SMEs is compensated by an increase 
in the size of the form factors. This unexpected result means that the obtained CPU 
performance — at least for this process — does not depend on sophisticated and process- 
dependent optimisations. 

2.2 Real corrections 

The calculation of the qq channel has been described in Ref. |TT]. The evaluation of 
the g'g channels is done in the same way. In the following we sketch the calculation for 
the gg channel. We have again performed two independent calculations of all building 
blocks. 

In both evaluations of the real corrections the amplitudes are calculated as helicity 
matrix elements which have been generated with Madgraph 4.1.33 [67j. While the 
amplitudes for qq — )• ttbbg and Vs ~^ ttbbVhave been checked with the Weyl-van der 
Waerden spinor formalism of Ref. [BS], those for gg — )■ ttbbg have been verified with 
an implementation of off-shell recursion relations p^HTTj. The singularities for soft and 
collinear gluon emission are isolated via dipole subtraction [501453] for NLO QCD calcu- 
lations using the formulation [53] for massive quarks. Soft and collinear singularities in 
the "endpoint part" of the subtraction function (the / operator of Refs. [50l|53]), i.e. the 
part of the subtraction terms that has to be combined with the virtual corrections, are 
regularized using the same regularization prescription (dimensional or with small quark 
masses) as the corresponding virtual corrections. No regularization is needed in the sub- 
traction terms for the real corrections. For both the qq and gg channels 30 different dipole 
subtraction terms need to be included while each gg channel requires only 10, since we 
demand b quarks with finite transverse momentum in the final state. After combining 
virtual and real corrections, singularities connected to collinear configurations in the final 
state cancel for "collinear-safe" observables after applying a jet algorithm. Singularities 
connected to collinear initial-state splittings are removed via MS QCD factorization by 
PDF redefinitions. In both evaluations the phase-space integration is performed with 
multi-channel Monte Carlo generators [51] and adaptive weight optimisation similar to 
the one implemented in LusiFER [56j. 

Version 1 of the real corrections employs the MadDipole implementation of dipole 
subtraction [72]. The phase-space integration, implemented in C-|--|-, is based on 
RacoonWW, but the phase-space mappings are built up in a more generic way very 

® Products of four Dirac matrices can occur only inside the massive top-quark chain. In this case 
one can introduce a fifth Dirac matrix by rewriting the massive spinor as u{p) — ^u{j))/m^ and then use 
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c^/clo 


^ events 
(after cuts) 


(Acr)stat/o" 


runtime 


time /event 


tree level 


86% 


5.3 X 10^ 


0.4 X 10-3 


38 min 


0.4 ms 


virtual 


-11% 


0.26 X 10*^ 


0.6 X 10-3 


13 h 


180 ms 


real + dipoles 


49% 


10 X 10^ 


3 X 10-3 


40 h 


14 ms 


total 


124% 




4 X 10-3 


53 h 





Table 1: Statistics and speed of various parts of the calculation based on 2 x 10 events 
before applying cuts generated on a 3 GHz Intel Xeon processor using the pgf77 compiler 
with standard options. 

similar to the approach of LusiFER [56]. For each of the 341 bremsstrahlung Feynman 
diagrams a corresponding channel is taken into account in the Monte Carlo integration. 

In version 2 all dipole subtraction terms have been implemented directly into the 
Monte Carlo generator. The Monte Carlo generator is a further development of the 
one used in COFFER77 [73] and for the calculation of the NLO corrections to pp — t- 
Hjj + X [7^ . In addition to the Monte Carlo channels for the bremsstrahlung diagrams 
30 X 36 = 1080 channels are used to map the dipole subtraction terms in the gg channel. 
These additional channels lead to some improvement in the convergence of the Monte 
Carlo integration. 

All real bremsstrahlung amplitudes of Madgraph have been checked against indepen- 
dent calculations for several phase-space points. The cancellation between real matrix el- 
ements and dipole subtraction terms has been verified numerically in all soft and collinear 
regions. The individual dipole subtraction terms, the subtracted real matrix elements, and 
the integrated subtraction terms (P and K terms of Refs. [SOIES]) have been compared 
point-wise between the two independent calculations. The agreement was generally at 
the level of 13 digits. The integrated LO cross section has been verified with SHERPA [75] 
at the level of the integration errors of 0.2%. From the point of view of Monte Carlo 
integration the most complicated and time-consuming part is the integration of the real 
corrections in the gg channel. For the complete NLO cross section we found agreement 
between the two versions of our code within l-3(j, where la corresponds to 0.1-0.2%. 
The results for the distributions coincide within 1-3(T for each bin. 

In order to give an impression on the statistics and the required CPU time we show in 
Table [T]some numbers for 2 x 10'' generated phase-space points before cuts. This yields an 
accuracy for the NLO cross section of about 0.5%. The contributions of the Vs channel 
were calculated for every 4th event and those of the virtual corrections in the gg channel 
and the contributions of the qq channel for every 20th event. The bulk of the runtime is 
taken by the gg channel. For the virtual corrections the CPU time is dominated by the 
gg channel and amounts to 180 ms per event. The virtual correction in the qq channel are 
by a factor 20 faster. In order to produce the plots for the scale variations we generated 
2 X 10'' phase-space points for the NLO predictions and 2 x 10*^ for the LO cross section. 
To generate the distributions we used about 20 times more events for the NLO results. 
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3 Predictions for the LHC 



The thorough description of a complex signature like ttbb involves numerous possible 
observables. Here we investigate distributions and cuts that are relevant for the search 
of ttH(H — )• bb) at the LHC [5l[7], where ttbb contributes to the irreducible background. 
In our previous work [TT] we found a K factor of about 1.8 for the integrated ttbb cross 
section at the LHC. This unexpectedly large NLO effect raises two important issues that 
we address in the present analysis. Firstly, we discuss the relation between the large K 
factor and the scale choice. This leads us to a new and more appropriate scale choice, 
which improves the convergence of the perturbative expansion. Secondly, we consider 
possible strategies to reduce the ttbb cross section in order to facilitate the extraction 
of the ttH signal. In particular, we study the influence of a jet veto on the NLO cross 
section and its perturbative stability. We also explore the kinematic region of highly- 
boosted bb pairs, which helps to separate the Higgs signal from its QCD background, as 
flrst suggested for associated WH and ZH production jSj and, very recently, also for ttH 
production [9]. 

Let us remind that top-quark decays are not included in our calculation. In practice we 
treat top quarks in a completely inclusive way, and we restrict our analysis to the kinematic 
properties of those two b quarks that do not result from top decays. From the experimental 
view-point, the presented distributions correspond to the unrealistic situation of perfect 
top-quark reconstruction. The detailed description of top-decay products and the related 
issue of b-quark combinatorics are left for future studies. 

3.1 Input parameters, jet definition, and cuts 

We study the process pp — ttbb -|- X at = 14 TeV. For the top-quark mass, 
renormalized in the on-shell scheme, we take the numerical value = 172.6 GeV [76] . 
All other QCD partons, including b quarks, are treated as massless particles. CoUinear 
flnal-state conflgurations, which give rise to singularities, are recombined into IR-safe 
jets using a fcx-algorithm [77]. Speciflcally, we adopt the fcx-algorithm of Ref. [78] and 
recombine all flnal-state b quarks and gluons with pseudorapidity \ri\ < 5 into jets with 
separation a/A</)2 + A?/^ > D = 0.4 in the rapidity-azimuthal-angle plane. Requiring two 
b jets, this also avoids collinear singularities resulting from massless g — )■ bb splittingsjll 

After recombination, we impose the following cuts on the transverse momenta and 
rapidities of the b jets: 

PT,h > PT,b,cut = 20 GeV, |?/b| < 2/b,cut = 2.5. (3.1) 

This choice is dictated by the detector geometry and the search for a ttH(H — )■ bb) 
signal at the LHC [SIE]. The outgoing (anti)top quarks are neither affected by the jet 
algorithm nor by phase-space cuts. For what concerns b quarks, the jet algorithm and the 
requirement of having two b jets with pT,b > 20 GeV sets an effective lower limit on the bb 
invariant mass of roughly 10 GeV. But the m^b-range that is relevant for the Higgs-boson 
search is actually much higher. In this kinematic region, with pT,b ^ f^h and m^b ^ ^b, 

^ Note that, as compared to our previous analysis [lOlIll) . we have reduced the jet- algorithm parameter 
from D = 0.8 to D = 0.4. This is particularly important for highly boosted b-quark pairs with TOj^g ~ Mh, 
since D = 0.8 would lead to their recombination into a single jet and, consequently, to their rejection. 
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we expect that the rrib = approximation works fairly well. To asses its precision we 
compared the LO cross section for = and my, = 4.2 GeV using SHERPA [75]. For the 
integrated cross section, which is dominated by m^b values well below 100 GeV, we found 
that the finite-rrib effect is smaller than 3%. 

We consistently use the CTEQ6 [79] set of PDFs, i.e. we take CTEQ6L1 PDFs with 
a one- loop running as in LO and CTEQ6M PDFs with a two-loop running as in NLO, 
but neglect the suppressed contributions from b quarks in the initial state. The number 
of active flavours is A^p = 5, and the respective QCD parameters are A^^ = 165 MeV and 
y^MS _ 226 MeV. In the renormalization of the strong coupling constant the top-quark 
loop in the gluon self-energy is subtracted at zero momentum. In this scheme, the running 
of as is generated solely by the contributions of the light-quark and gluon loops. 

3.2 Renormalization and factorization scales 

The perturbative expansion of the pp — )■ ttbb cross section starts with the fourth power 
of as. Consequently the LO predictions, and also the K factor, are extremely sensitive to 
variations of the renormalization scale. In Ref. [1] it was pointed out that the LO ttbb 
cross section can vary by up to a factor four if the QCD scale is identified with different 
kinematic parameters. In all recent ATLAS studies of ttH(H — )■ bb) [2i|4|l5] the signal 
and its ttbb background were simulated by setting the renormalization and factorization 
scales equal to half the threshold energy, ii^thr = '^^t + ^hh- For pp — )■ ttH-l-X, this choice 
was well supported by the existing NLO analysis [T3l[Tl]. But, in the absence of NLO 
predictions for ttbb, the choice of the same scale for signal and background was motivated 
solely by the assumption that the two processes have similar kinematics. However, in 
Ref. [11] we found that, if both processes are evaluated at /ir = /^f = -£'thr/2, pp — )■ ttbb 
receives much larger NLO corrections {K ~ 1.8) than pp — )■ ttH {K ~ 1.2). This is 
mainly due to the fact that the scale -Ethr/2 does not provide an adequate description of 
the QCD dynamics that governs ttbb production. 

The main difference between pp — > ttH(H — )■ bb) and its irreducible QCD background 
is that the former process involves only two powers of as at LO. Moreover, the part 
of the signal process that is mediated by strong interactions does not involve any scale 
significantly smaller than i?thr/2. In contrast, pp — t- ttbb is entirely driven by QCD and is 
proportional to at LO. In the m^b — ^ limit, the dominant ttbb production mechanism 
is pp — )■ ttg{g — )■ bb), where a gluon with small virtuality plays a role analogous to the 
intermediate Higgs boson in the signal process. In this regime, the factorization of ttg 
production and g — t- bb splitting provides two well-defined and widely separated scales: 
£'thr/2 and mbS, respectively. This simple picture is, however, not applicable to the 
kinematic region of interest, where mbb ^ 100 GeV. Here, pp — )■ ttbb involves various 
other mechanisms. For instance, the radiation of one or both b quarks off initial-state 
gluons can play an important role due to collinear enhancements. In order to find an 
optimal QCD scale, we have tried to identify a dominant production mechanism. To 
this end, we have inspected the relative weights of the channels corresponding to various 
Feynman- diagram topologies in our adaptive Monte Carlo generator. However, we found 
that none of these channels is strongly enhanced with respect to the others. Similarly we 
were not able to reproduce the bulk of the cross section in terms of effective approximations 
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Setup 


^bb,cut PT,bb,cut Pjet,veto 


PT,b,cut |/b,cut 


clo ctnlo K 


I 

II 
III 
IV 


100 

200 

100 - 100 
100 


20 2.5 
20 2.5 
20 2.5 
50 2.5 


786.3(2)+^^^; 978(3) 1.24 
451.8(2)^^^^; 592(4) 1.31 

786.1(6) +2% 700(3) 0.89 
419.4(1) 526(2) 1.25 



Table 2: Cut parameters (in GeV), integrated LO and NLO cross section (in fb) with 
statistical errors and scale variations by factors 2 up and down as well as K factors for 
the four different setups. 



based on collinear g — bb splittings of the incoming gluons. This suggests that pp — > ttbb 
is a genuinely multi-scale and multi-channel reaction. 

We thus decided to adopt a pragmatic scale choice, based on the kinematic properties 
of the ttbb final state. While rrit sets a clear scale for the couplings to the top quarks, 
the inspection of differential distributions reveals that the cross section is saturated by b 
quarks with px.b ^ f^t (see Figs. |8]and|9]). In order to account for these different scales 
we have adopted a dynamical QCD scale corresponding to their geometric average, 

/^o = mt^yp^,hPT,h■ (3-2) 

Our LO and NLO predictions are obtained by varying the renormalization (/xr) and 
factorization (fip) scales around the central value fl3.2p . 

yWR = ^RfJ'O, fJ'F = ^F/io- (3.3) 

In the following sections we investigate the dependence of the LO and NLO integrated 
cross section with respect to uniform (^p = ^r) and antipodal (,^f = ■Cr^) scale variations 
in the range 1/8 < C,f,C,b. < 8. We find that uniform variations have a larger impact on 
the cross section as compared to antipodal variations. For all distributions we provide 
LO and NLO predictions with uncertainty bands corresponding to factor-two uniform 
scale variations. More precisely, all observables are evaluated at three different scales: 
^F = ^R = 0.5,1,2. As we will see from the reduction of the K factor and the scale 
uncertainties, the scale choice ( 13. 2 p clearly improves the convergence of the perturbative 
expansion as compared to Ref. [11]. 

3.3 Additional cuts 

Besides the standard cuts (13. ip . we have imposed the following kinematic restrictions 
to the bb system and the extra jet that is radiated at NLO: 

m^b > "^bb.cut) PT,bb > PT,bb,cut5 PTJet < Pjet,veto- (3.4) 

In order to investigate the individual effect of these extra cuts and correlations with 
other observables, we have generated differential distributions in four different setups 
described in Table [2l The setups I-III implement the standard cuts ( 13. ip and explore the 
individual impact of the extra cuts (13. 4p . Setup IV is a variant of setup I, where the cut 
( 13. ip on the b-jet is increased from 20 GeV to 50 GeV. 
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Figure 4: Scale dependence of the LO and NLO pp — )■ ttbb + X cross section at 
^/s = 14TeV in setup L The left and the right plots describe uniform (^r = = 
and antipodal (^r = = ^) scale variations, respectively. 



3.4 Setup I 

In this setup we impose the cut m^s > 100 GeV. This removes a large part of the 
cross section and selects the kinematic region of interest for the ttH(H — > bb) signal. 

Scale dependence 

The LO and NLO integrated cross sections and their dependence with respect to 
uniform (left plot) and antipodal (right plot) scale variations are displayed in Figure HI 
At the central scale we obtain o"lo = 786.3(2) fb and o"nlo = 978(3) fb, where the numbers 
in parentheses are the errors of the Monte Carlo integration for 2 x 10® LO events and 
2 X 10'' NLO events before applying cuts. These predictions are not directly comparable 
to those of Ref. [11], where we did not apply any cut to m^^. Still we can compare the 
K factors, which are rather insensitive to m^^. We observe that the new scale choice 
(13. 2p reduces the NLO corrections from K ~ 1.77 [llj to ~ 1.24. We note that, 
in spite of the smaller K factor, the new scale choice yields larger LO and NLO cross 
sections as compared to the scale -Ethr/2 used in Ref. [H]. This can be easily seen from 
Figure m where the new and the previous central scales correspond to /iR//io = 1 and 
Mr/Mo = -Ethr/(2yUo) ^ 1, respectivelylf] In addition to the K factor, the new scale choice 
reduces also the NLO scale uncertainty. Varying the scale up or down by a factor 2 changes 
the LO and NLO cross section by 78% in LO and by 21% in NLO. The improvement with 

^ Using setup-I cuts and the ATLAS scale choice ^r^f = £^thr/2 we obtain crLo(-Ethi/2) = 448.7(1) fb 
and (TNLo(£'thi72) — 751(2) fb. The increase of the cross section due to the combined effect of the new 
scale choice and the NLO correction factor is thus crNLo(/^o)/o'Lo(£'thi72) ~ 2.18 while the K factor for 
setup-I cuts and the ATLAS scale choice results in if ~ 1.67. 
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Pjet,veto [GeV] Pjet,veto [GeV] 

Figure 5: Dependence of the pp — )■ ttbb + X cross section on a jet veto (pxjet < Pjet,veto) 
in setup I: absolute LO and NLO predictions (left) and NLO K factor (right). The 
uncertainty bands correspond to factor-two scale variations. 

respect to Ref. [H], where we had a 33% NLO uncertainty, is evident also from the shape 
of the NLO curves in Figure HI Now we observe a stable point in the vicinity of the central 
scale. 

Jet veto 

As anticipated in Ref. [11], a jet veto can significantly reduce the large cross section 
of the ttbb background. This could facilitate the extraction of the ttH signal at the 
LHC. In Figure [H the integrated ttbb cross section is plotted versus the upper bound, 
Pjet.veto, imposed to the jet transverse momentum. Here, as well as in the following figures, 
the left plot shows the absolute LO (blue) and NLO (red) predictions. The curves and 
their uncertainty bands represent factor-two (uniform) scale variations around the central 
value (13.21) . The right plot displays the LO and NLO bands normalized to the central 
value of the LO prediction. There the blue band illustrates the relative uncertainty of 
the LO cross section, i.e. crLo(0/^Lo(^ = 1), the red curve corresponds to the K factor, 
K = crNLo(^ = l)/c"Lo(^ = 1), and the red band shows the variation of the K factor 
when varying the scales in the NLO cross section but keeping them fixed in the LO cross 
section, = ctnlo (0 /"^lo (^ = !)• Figure [5] the red (NLO) curve tends to saturate 
the upper bound of its uncertainty band, a feature that can be observed in various other 
distributions. This is consistent with the shape of the NLO curve in Figure HJ which 
develops a maximum in the vicinity of the central scale. 

The NLO curve in Figure [5] shows that a sizable reduction of the cross section requires 
a jet veto well below 200 GeV. For pjet.vcto = 150, 100, and 50 GeV, the K factor is reduced 
to 1.03, 0.89, and 0.54, respectively. However, there is a trade-off between suppressing 
the NLO cross section and increasing its perturbative uncertainty. The jet veto tends to 
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destroy the cancellation between IR logarithms of virtual and real origin and its effect 
grows as — ln^(i?thr/Pjet,veto) when pjet.veto becomes small. Since they are accompanied 
by an coefficient, these logarithms can give rise to huge scale uncertainties already for 
moderate values of Pjet,veto- 

This is reflected by the dramatic amplification of the NLO uncertainty band in 
Figure [51 Its lower bound enters the pathologic regime of negative cross sections around 
Pjet,veto = 50 GeV. Here, besides the NLO cross section itself, also its uncertainty es- 
timate becomes completely unreliable. The region of small Pjet,veto would require a re- 
summation. This would stabilize the perturbative calculation and compensate for its 
divergent behaviour. As a result, the unphysical suppression of the NLO cross section 
for Pjet.veto ^ 100 GeV would be washed out. If we restrict ourselves to the fixed-order 
NLO result, the plot tells us that jet- veto values around 100 GeV provide a good compro- 
mise: the reduction of the K factor is already significant {K ~ 0.89) and the NLO scale 
uncertainty (19%) is at the same level as for the total cross section (21%). 

Invariant-mass distributions 

The invariant-mass distribution of the bb pair, shown in Figure |6l constitutes a key 
observable for the search of the ttH signal. Because of limited resolution and b-quark 
combinatorial problems, the Higgs boson would appear as a relatively broad and small 
peak on top of this distribution. The subtraction of the ttbb background requires an 
accurate determination of its normalization, possibly by direct measurement in a signal- 
free region, and a precise theoretical description of its shape. 

In Figure |6] we observe that the NLO predictions perfectly fit within the LO band 
and significantly reduce the QCD uncertainty over the entire invariant-mass range. The 
numerical impact of the corrections is moderate and almost constant (1.21 < K < 1.27). 
This favourable behaviour is ensured by the dynamical scale choice f l3.2p . At high invariant 
masses the upper bound of the NLO uncertainty band slightly decreases and approaches 
the central NLO prediction. The same trend appears in the high-energy tail of other 
distributions. 

Figure [7] displays the dependence of the cross section with respect to a cut on the 
total invariant mass of the ttbb system {m^i^^ > m^ut.ttbb)- Since it corresponds to the 
invariant mass of the W"^W~bbbb final state, this quantity is independent of the b-jet 
combinatorics. It can thus be measured with better resolution as compared to observables 
that involve only a particular subset of the four b jets. This property might be exploited 
in order to improve the signal-to-background ratio. Apart from a slight increase in the 
high invariant-mass tail, the NLO corrections behave similarly as for the m^\^ distribution. 

Transverse-momentum distributions 

The transverse-momentum distributions of the individual b jets, ordered according to 
their p^, are presented in Figure |8] (harder b jet) and Figure [9] (softer b jet). While the 
softer b jet tends to saturate the lower bound of 20 GeV (13. ip . the harder behaves rather 
differently. Its distribution has a maximum around 100 GeV and a tail that extends up to 
fairly high transverse momenta. These shapes suggest that one of the two quarks is often 
emitted from initial-state gluons, while the other one participates to the hard scattering. 
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Figure 6: Invariant-mass distribution of the bb pair in setup I: absolute LO and NLO 
predictions (left) and NLO K factor (right). The uncertainty bands correspond to factor- 
two scale variations. 
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Figure 7: Dependence of the cross section with respect to a cut on the ttbb invariant 
mass (m^tbh > ''^ttbb,cut) ill setup L absolute LO and NLO predictions (left) and NLO K 
factor (right). The uncertainty bands correspond to factor-two scale variations. 
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Figure 8: Transverse-momentum distribution of the harder b jet in setup I: absolute LO 
and NLO predictions (left) and NLO K factor (right). The uncertainty bands correspond 
to factor-two scale variations. 
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Figure 9: Transverse-momentum distribution of the softer b jet in setup I: absolute LO 
and NLO predictions (left) and NLO K factor (right). The uncertainty bands correspond 
to factor-two scale variations. 
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Figure 10: Transverse-momentum distribution of the bb system in setup I: absolute LO 
and NLO predictions (left) and NLO K factor (right). The uncertainty bands correspond 
to factor-two scale variations. 



In contrast, none of the b quarks resulting from ttH originates from initial-state radiation. 
This feature, which renders the cross section quite sensitive to pT,b, might be exploited 
to improve the separation of the ttH signal. 

The dynamical scale introduced in (13.21) accounts for the different kinematics of the 
two b jets and the extension of their transverse momenta over a wide pt range. The 
goodness of this choice is confirmed by the stability of the K factor over the entire px- 
spectrum. The transverse-momentum distribution of the bb pair is shown in Figure [TOl 
Its shape resembles fairly closely that of the harder b-jet distribution. Also the K factor 
and the scale uncertainties behave similarly. 



Rapidity and azimuthal distributions 

The rapidities of the individual b jets, ordered in px, are shown in Figure ITT] 
(harder b jet) and Figure fT2] (softer b jet). Both b jets tend to populate the central 
region. But this feature is much more pronounced for the harder b jet, while the softer 
one has a significant probability to be emitted also in the forward and backward direc- 
tions. The rapidity distribution of the bb system (not plotted) resembles that of the 
harder b jet, and the rapidity-separation distribution (Figure [13]) does not suggest any 
strong correlation between the two b jets. All rapidity distributions receive moderate and 
almost constant NLO corrections. 

The rapidity-azimuthal-angle separation, Ai?bb = \J {Uh — UhY + (0b — 0b)^! of the 
b jets is displayed in Figure [TH The shape of this distribution is determined by three 
kinematic constraints: the rapidity cut (13. ip is responsible for the suppression at high 
Ai^bs; the sharp lower bound at Ai^^B = 0.4 results from the jet algorithm; and the 
invariant-mass cut m^b > 100 GeV keeps the two b jets at intermediate Ai?-separations. 
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Figure 11: Rapidity distribution of the harder b jet in setup I: absolute LO and NLO 
predictions (left) and NLO K factor (right). The uncertainty bands correspond to factor- 
two scale variations. 
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Figure 12: Rapidity distribution of the softer b jet in setup I: absolute LO and NLO 
predictions (left) and NLO K factor (right). The uncertainty bands correspond to factor- 
two scale variations. 
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Figure 13: Rapidity separation of the two b jets in setup I: absolute LO and NLO predic- 
tions (left) and NLO K factor (right). The uncertainty bands correspond to factor-two 
scale variations. 
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Figure 14: Rapidity-azimuthal-angle separation of the two b jets in setup I: absolute LO 
and NLO predictions (left) and NLO K factor (right). The uncertainty bands correspond 
to factor-two scale variations. 
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Figure 15: Azimuthal orientation of the b jets in the plane perpendicular to the bb 
(see text) in setup I: absolute LO and NLO predictions (left) and NLO K factor 
The uncertainty bands correspond to factor-two scale variations. 



system 
(right). 



COS0bb = 



The NLO corrections induce a 30-40% distortion of the shape of this distribution in the 
region 0.4 < Ai?bb ~ 1- This effect can be attributed to the recombination of b quarks 
and non-b partons, which can turn b-quark pairs with AR < 0.4 into b-jet pairs with 
AR > 0.4. 

Finally, in Figure [T5| we plot the angular variable 0bb5 which describes the azimuthal 
orientation of the b jets. This observable is defined as the opening angle between two 
planes: the "production" plane spanned by the beam axis and the total momentum of 
the bb system, and the "decay" plane, which contains the momenta of the individual b 
jets, defined by 

[Pbcam X (Pb + Pb)] ■ (Pb X Pb) ,^ 
IPbeam X (pb + Pb)||Pb X Pb| ' 

Equivalently 0bb represents the azimuthal orientation of the b jets with respect to the 
beam direction in the plane perpendicular to the bb momentum. In the case where the 
bb pair results from the decay of an intermediate particle, like the Higgs boson in ttH 
production, the spin of the latter can be determined from the 0bb distribution. Since 
the Higgs boson has spin 0, the 0bb distribution is expected to be isotropic, i.e. 0bb- 
independent. As we see from Figure fT5| the NLO corrections have a negligible influence 
on the shape of this observable for the ttbb background. 

3.5 Setup II 

As discussed in the introduction, the selection of ttbb signatures with highly boosted 
b-quark pairs may help to separate the ttH(H — )• bb) signal from its backgrounds. This 
motivates us to study the irreducible ttbb background in this particular phase-space 
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Figure 16: Scale dependence of the LO and NLO pp — )■ ttbb + X cross section at 
^/s = 14TeV in setup IL The left and the right plots describe uniform = = 
and antipodal (,^r = = ^) scale variations, respectively. 



region. Specifically, in setup II, we select highly boosted bb pairs with px.bb > 200 GeV, 
as proposed in Ref. [9]. In contrast to setup I, here we do not impose any cut on the bb 
invariant mass. Nevertheless the cuts on Px.bB) PT,h, and Px.S^ together with the bound 
ARy^^ > 0.4 resulting from the jet algorithm, impose an effective lower bound 

^ Ai?bbPT,bb - z) > Ai?bbPT,b,cut r'^'^^'^"' - 1 = 24 GeV, (3.6) 

where z and 1 — z are the (transverse) momentum fractions of the two b jets, and the 
first equation holds for small Ai^^g. 



Scale dependence 

The scale dependence of the LO and NLO integrated cross sections is shown in 
Figure [TBI At the central scale we obtain ctlo = 451.8(2) fb and (Jnlo = 592(4) fb. This 
corresponds to an NLO correction factor K = 1.31. The absolute NLO cross section is 
reduced by about 40% as compared to setup I. The shape of the scale-dependence curves 
is quite similar as in Figure H] and indicates good convergence and stability of the pertur- 
bative expansion. The shifts induced by factor-two variations of the QCD scales amount 
to 79% in LO and 22% in NLO. 

Investigating the sensitivity of the NLO cross section to a jet veto we found similar 
results as in setup I. For a jet veto of 100 GeV the K factor and the NLO uncertainty 
amount to 0.84 and 33%, respectively. 
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Figure 17: Invariant-mass distribution of the bb pair in setup II: absolute LO and NLO 
predictions (left) and NLO K factor (right). The uncertainty bands correspond to factor- 
two scale variations. 



Invariant-mass and transverse-momentum distributions 

The bb invariant-mass distribution is displayed in Figure [T71 Its behaviour in the 
region m^b ^ 50 GeV reflects the effective lower bound (13.61) . We find that the NLO 
corrections induce an appreciable shape distortion of about 20%, in particular near the 
physically interesting region of m^b ~ 100 GeV. Such effect tends to mimic a Higgs signal 
and should be carefully taken into account in the ttII(II — )■ bb) analysis. 

The transverse-momentum distributions of the harder and softer b jets are presented 
in Figure [18] and Figure [191 respectively. As a consequence of the cut imposed on the 
transverse momentum of the b pair, the harder b jet is pushed to much higher values 
as compared to setup I. The maximum of its distribution is located around 200 GeV. 
In contrast, the softer b jet is much less sensitive to the Px.bB cut and is predominantly 
produced in the region 20 GeV < < 100 GeV. This different kinematic behaviour of 
the two b jets might be exploited to separate the ttbb background from the ttH signal, 
where both b jets are produced by the Higgs boson and should thus have more similar 
PT-values. The NLO corrections to both pT,b distributions feature a slight transverse- 
momentum dependence, with 10% variations of the K factor within the plotted range. 

Rapidity and azimuthal distributions 

The rapidities of the harder and softer b jets are shown in Figures [20] and [211 respec- 
tively. Due to the high of the b-pair system, both b jets tend to be more central as 
compared to setup I. While the K factor is almost insensitive to the rapidity of the soft 
b jet, the NLO corrections have a non-negligible influence on the shape of the hard-b-jet 
distribution for \y\,\ > 1. The rapidity distribution of the bb system (not plotted) behaves 
similarly to the one of the harder b jet. 
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Figure 18: Transverse-momentum distribution of the harder b jet in setup II: absolute LO 
and NLO predictions (left) and NLO K factor (right). The uncertainty bands correspond 
to factor-two scale variations. 
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predictions (left) and NLO K factor (right). The uncertainty bands correspond to factor- 
two scale variations. 
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Figure 22: Rapidity separation of the two b jets in setup II: absolute LO and NLO 
predictions (left) and NLO K factor (right). The uncertainty bands correspond to factor- 
two scale variations. 



The rapidity-separation distribution (Figure [22]) is strongly peaked at small Ay^b and 
the NLO corrections have an appreciable influence on its shape. In the region Ay-^Y, < 2 
the K factor varies between 1.17 and 1.38. 

Finally, in Figure [23] we show the rapidity-azimuthal-angle separation of the b jets, 
which is strongly peaked at small AR^^,. Also in this distribution we observe a significant 
shape distortion. In the region AR^^, < 1, which corresponds to invariant masses of the 
order of 100 GeV, the K factor increases by about 20%. This NLO effect might have an 
important impact on the measurement of ttH in the highly-boosted Higgs regime. 

3.6 Setup III 

In order to explore the effect of a jet veto and its possible correlation with other 
observables, we have generated events with m^s > 100 GeV and PTjet < 100 GeV. The 
LO and NLO cross sections and their scale dependence are shown in Figure [24] At the 
central scale we obtain ctlo = 786.1(6) fb and ctnlo = 700(3) fb, corresponding to a 
correction factor = 0.89. It is evident from Figure [2^ that the central scale is very 
close to a stable point. This demonstrates that a jet-veto value of 100 GeV is sufficiently 
large to avoid perturbative instabilities. Varying the QCD scales up and down by a factor 
two shifts the NLO cross section by only 0.4% and —19%, respectively. 

Inspecting various kinematic distributions we find that the NLO corrections have a 
much bigger impact on shapes as compared to setup I. In particular, we observe that 
the suppression effect resulting from the jet veto is rather sensitive to the transverse 
momentum of the b jets. For instance, the NLO correction to the pt distribution of the 
harder b jet varies by about 15% in the range 20 GeV < PT,bi < 200 GeV. 
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NLO predictions (left) and NLO K factor (right). The uncertainty bands correspond to 
factor-two scale variations. 
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Figure 24: Scale dependence of the LO and NLO pp ttbb + X cross section at 
y/s = 14TeV in setup III. The left and the right plots describe uniform (^r = (,f = 
and antipodal (^r = = ^) scale variations, respectively. 
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Figure 25: Scale dependence of the LO and NLO pp — )■ ttbb + X cross section at 
= 14TeV in setup IV. The left and the right plots describe uniform = = 
and antipodal (,^r = = ^) scale variations, respectively. 



3.7 Setup IV 

As observed in Figure [HI the ttbb cross section is dominated by relatively soft b jets, 
which tend to saturate the pT,b > 20 GeV cut. It is thus interesting to investigate the 
influence of this cut on the ttbb background and, in particular, on the NLO corrections. 
To this end, we have studied a variation of setup I where the pT,b cut is increased from 
20 to 50 GeV. The LO and NLO cross sections and their scale dependence are shown 
in Figure [251 At the central scale we obtain ctlo = 419.4(1) fb and ctnlo = 526(2) fb, 
corresponding to a correction factor K = 1.25. The higher pT,b cut reduces the NLO 
cross section by 46% as compared to setup I. The LO and the NLO scale dependence 
remain very similar as in setup I (cf. Figure H]): the uncertainty corresponding to factor- 
two scale variations amounts to 77% in LO and 21% in NLO. 

Inspecting various kinematic distributions we find, as in setup I, that the NLO correc- 
tions have little impact on the shapes. In general the kinematic dependence of the NLO 
correction factor is below 10%. The largest shape distortion is observed in the AR^^^ 
distribution: similarly as in setup I (cf. Figure [T4|) we find an increase of the K factor in 
the region 0.4 < Ai^^b ~ 1 by about 30%. 



4 Conclusions 

The direct production of ttbb final states represents the major background to the 
production of Higgs bosons in association with top-antitop-quark pairs at the LHC, 
pp — 7- ttH, where the Higgs boson decays into bb pairs. This process can lead to a 
direct measurement of the top-quark Yukawa coupling. Apart from improvements in the 
experimental analysis, a successful exploitation of this demanding channel requires predic- 
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tions for ttbb production at the next-to-leading-order level in QCD, maybe even further 
improvements. 

Extending our first results on the total cross section for pp — > ttbb pubhshed earlier, 
wc have presented a detailed study of integrated and differential cross sections at the 
LHC, discussing in particular a dynamical scale choice, the influence of various cuts on 
the outgoing b quarks, and the impact of a veto on an additional hard jet. We observe 
that the traditional choice of a constant scale determined by the energy threshold for 
the process underestimates the ttbb background by a factor of two, while an appropriate 
dynamical scale, which is tied to the transverse momenta of the b quarks, stabilizes the 
perturbative predictions much better. Moreover, the K factor is reduced from 1.8 to 1.2. 
Using the new scale choice, the corrections have little impact on the shapes of distributions 
if standard cuts are applied. Strengthening the cut on the transverse momentum of the 
bottom quarks has no big influence on the K factor and the effect of the NLO corrections 
on the shape of the distributions. On the other hand, imposing a jet veto of 100 GeV 
reduces the K factor to 0.9 and enhances the impact of the NLO corrections on the 
shapes of distributions. In the regime of highly boosted Higgs bosons, which offers better 
perspectives to observe the ttH signal, we find significant distortions of the kinematic 
distributions, while the K factor is 1.3. 

Our calculation builds on the Feynman-diagrammatic approach, i.e. on an algorithmic 
reduction of each Feynman diagram to a canonical standard form, which is automatically 
processed to FORTRAN output, and on a numerical reduction of tensor loop integrals to 
an appropriate set of scalar master integrals. A key feature of the diagram-by-diagram 
approach is that colour sums can be preformed very efficiently. Helicity summation is 
simplified by introducing a basis of 0(1000) basic structures. The reduction to these 
structures can be performed in a process-independent way. The numerical tensor-integral 
reduction employs dedicated methods that have been developed to treat the numerically 
delicate phase-space regions where small Gram determinants appear in denominators 
during the traditional tensor reduction. Real corrections are integrated using well-known 
dipole subtraction methods. Wc find that our Feynman-diagrammatic approach provides 
very high numerical stability and CPU efficiency, a result that is very encouraging in 
view of future challenging next-to-leading-order calculations for important multiparticle 
processes. 
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Appendix 



A Some details on the reduction of standard matrix elements 

In Section [2?T] we have briefly described our prodecure to reduce the numerous helicity 
structures for gg — )■ ttbb to a standard form. We proceed in two steps, the flrst step 
employing only identities that hold in arbitrary D ^ 4 space-time dimensions, followed 
by the second step (after cancelling UV divergences) that builds on four- dimensional 
identities. As also anticipated in Section |2TT| we support two different versions of step 2. 
The flrst variant only uses the reduction (12.211) of products of flve Dirac matrices to 
products of three or only single matrices, leading to 970 SMEs. In this approach no chiral 
projectors or 75 factors are introduced. The second variant employs more four- dimensional 
identities, such as (I2.20p . which are derived from Chisholm's identity in Ref. [I0l[66], and 
leads to 502 SMEs, which however involve the 75 matrix. 

In the following we describe the second procedure for the crossed process 

g(Pi) g(P2) t(p3) t(p4) b(p5) Hpe) ^ 0, (A.l) 

where the incoming momenta of the corresponding incoming particles are given in paren- 
theses. The individual steps to reduce the structures 

[v{p3)TaUJaU{p4)][v{p5)Tb(^rU{pe)] = [Ta]s4[Tb\le (A.2) 

are as follows: 

1. First, we eliminate multiple contractions of Dirac matrices between the two Dirac 
chain^ by identities like 

7''7"7^w± ® -f^j^j^ = 4g''^Yi^± ® 7/.^^± + ^i'^uj± ® T^^^T' 

7''7"7''u;± ® -ful'^-ff, = %"'^7^w± ® l^i^T + ^1^^± ® 7"w±- (A.3) 

2. Next we shorten strings of flve or more Dirac matrices using (I2.2ip in the massless 
bottom-quark chain, \v(j)^Y\pj.ru{j)^\^ leaving three or only one Dirac matrix in this 
chain. There is at most one contraction of Dirac matrices between the two chains. 

3. Now we simplify structures that involve bottom-quark Dirac chains with three 
slashed vectors using the trick described in (3.40) and (3.41) of Ref. [66]. In more 
detail, this manipulation replaces the product of two Dirac chains of the form 
[(^^^ . . .]34[(^^/]5g by products in which either the top chain [. . .]34 has two Dirac 
matrices less or in which the bottom chain contains only one slashed vector. Us- 
ing this procedure recursively we obviously achieve that bottom chains with three 
slashed vectors are multiplied by top chains containing only a single Dirac matrix, 
i.e. they appear in the form [5^] 34 [5^^/] 56 • 

^In this context we recall that multiple contractions of Dirac matrices inside a single Dirac chain have 
already been eliminated in the first (ZJ-dimensional!) step of the algebraic reduction. 
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4. In a further step it is possible to eliminate all products of three Dirac matrices in the 
bottom chain, so that this massless Dirac chain contains exactly one Dirac matrix. 
If there is a contraction between the top and bottom chain, two Dirac matrices 
can easily be shifted from the bottom to the top chain using f l2.20p . If there is no 
contraction the product of chains looks like owing to the previous step. 

If the set {(i, e, /} contains at most one polarization vector Si {i = 1, 2) of the two 
incoming gluons, we can assume that {d, e, /} contains either the momentum or 
P4. If this is not the case, it can be achieved upon using momentum conservation 
and the Dirac equation in the bottom chain. The factors {i = 3, 4) in the bottom 
chain easily allow for shifting all but one Dirac matrix to the top chain by the 
manipulations. 
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where we used = {^j,7'^}/2 and the Dirac equation in the 34-chain and fl2.20l) . 
We recall that we consistently take momenta as incoming. The case that no factor 
1^3/4 appears in the bottom chain can only occur if the other two slashed vectors 
in the bottom chain belong to polarization vectors. Thus, we can assume that the 
vector a is a momentum out of {pi,P25P55P6} and that {rf, e, /} contains pi or p2- 
We can even omit p2 on either side if we eliminate it via momentum conservation, 
so that the considered bottom chains all contain a factor Two Dirac matrices 
can now be shifted from the bottom to the top chain using the relations 
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where we have used that the bottom mass is set to zero, i.e. = mg = 0. Note 
that we did not use = in the first relation, in order to indicate that the whole 
reduction procedure described here does not only apply to gluons but also to the 
case where gluons are replaced by massive vector bosons. 

5. After having reduced all bottom chains to contain only one Dirac matrix we apply 
f l2.2ip to the top chain recursively as far as possible. In this way the top chains 
contain only up to four Dirac matrices. 



6. Next we reduce products of the form 



34 LPJ56' 



which do not involve a contraction. 



After eliminating p2 via momentum conservation, the product contains one of 
the factors '^i'^^, '^I'^e, or "^^^^q, because at most two of the slashed vectors can be 
polarization vectors. The majority of such cases can be reduced with the following 
relations. 
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Here F is any string of Dirac matrices, and again ms = mg = was used. If a 7^ e, 
these relations interchange a and e (if a = e they are useless). Since e can only be 
Pi, P3, or p4, we can use these identities to shift ps or ^4 to the top chain, where these 
momenta are eliminated by the Dirac equation, so that we are left with the cases 
e = pi. For [F|^5|^6]34[|^i]56 this relation can be used to transfer pg to the bottom 
chain, again triggering the Dirac equation there. We are left with the two cases 
[F|^i|^5]34[|^i]5g and [F|^i|^g]34[|^i]5g, which are reduced according to 
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where the Dirac structure F and the vector a can be chosen arbitrarily. 



7. Now we reduce products of Dirac chains of the type [7^ 



34(7'^] 56; 



i.e. four Dirac 



matrices in the top chain with a contraction to the bottom chain. Among the vectors 
a, 6, c, there is at least one of the momenta pi, ps, pe- The momenta p5, pe can be 
shifted to the bottom chain via fl2.20p and subsequently eliminated with the Dirac 
equation. If only pi (but not P5, Pe) occurs in the top chain, the full Dirac structure 
can only be ['J^i'^i^i^2]M[l'^]56y which is kept as a standard structure. 

8. Now we repeat steps [6] and [7] for top chains involving three Dirac matrices, in order 
to eliminate such structures as far as possible, and subsequently again for top chains 
with two matrices. 

9. The remaining products of the form [|^i]34[|^j]56 are reduced as far as possible as 
described in (B.4) of Ref. [10]. 

10. Finally, it is convenient to replace all chirality projectors u± = {1 ± 75)/2 in terms 
of vector (no 75) and axial- vector (oc 75) structures in the two Dirac chains, because 
only the combinations vector®vector and axial-vector®axial-vector contribute ow- 
ing to the parity symmetry of the process in NLO QCD. 

B Benchmark numbers for the virtual corrections 



In order to facilitate a comparison to our calculation, in this appendix we provide 
explicit numbers on the squared LO amplitude and the corresponding virtual correction 
for a single non-exceptional phase-space point. The set of momenta for the partonic 
reaction gg — )■ ttbb is chosen as 



(500,0,0,500), 
(500,0,0,-500), 



p; 

= (327.5045589027869,107.1276753641986,-107.9290580423663,-233.1168284428635), 

P^ = (276.6425142763093,-107.4949148022111,153.8289259355409,-107.3397668261919), 

Pb = (233.9459027189062,82.55875671042013,-77.70592645955253,204.6375480757531), 

Ph 



(161.9070241019976,-82.19151727240762,31.80605856637796,135.8190471933023), 



fB.l) 



which defines the same phase-space point as already chosen for qq — )■ ttbb in 
The components are given in GeV and mt = 172.6 GeV. For the spin- and colour- 
squared LO amplitude we obtain 

IM^^l^gl = 5.437061775267626 ■ 10-^ GeV-^ 



Ref. [TO]. 

•averaged 



35 









5(1) 


5(0) 


loops 


version 1 
version 2 


-0.1484719139263099 
-0.1484719139260971 


0.0986990876957258 
0.0986990876958276 


0.4123948722195028 
0.4123948722188028 


I 


version 1 
version 2 


0.1484719139263437 
0.1484719139263439 


-0.0986990876957289 
-0.0986990876957291 


-0.2125642646365643 
-0.2125642646365644 


NLO 


version 1 
version 2 


0.0000000000000338 
0.0000000000002467 
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Table 3: Various contributions to the virtual NLO corrections to gg — > ttbb at the phase- 
space point (IB.ip . 



Madgraph 
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where we divided out the strong coupling constant g^. 

We express NLO contributions in the 2 — )■ 4 phase space as Laurent series in e 
(4-D)/2, 



\ fc=0 



(fc) , 



\M 



LO|2 



(B.3) 



where we factor out the LO term and the normalization factor 

(47r)T(l + e)r2(l-e) (4^ 



cr 



r(l-2e) 



r(l-e) 



+ 0{e 



IT 



6 



We split the result into the two parts, 

x(fe) _ 

"NLO ~ 



r(fc) 
"^loops 



(fc) 



(B.4) 



(B.5) 



which correspond to the contributions of renormalized loop diagrams (loops) and the I 
operator of the dipole subtraction function as defined in Ref . [53] . The numbers in Table [3] 
have been obtained in the 't Hooft-Feynman gauge using the 't Hooft-Veltman variant 
of dimensional regularization (four-dimensional external partons). We set the scale of 
dimensional regularization and the renormalization scale as /i = /iR = mt. The corre- 
sponding values of the strong coupling constant in the renormalization scheme described 
in Section [3] are 



as(^t)|LO = 0.1178730139006150, a,(mt)| 



NLO 



0.1076396017050965. 



(B.6) 



The agreement between our two independent versions of the virtual corrections is typically 
about 10 digits at regular phase-space points. 
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Another benchmark resuh for the virtual corrections to gg — )■ ttbb can be found in 
Ref. [H]. The quantity presented there and denoted as 'HELAC-IL' corresponds to the 
unrenormahzed virtual one-loop correction plus the mass-renormalization contribution 
(without wave-function and coupling constant renormalization) . We find agreement at 
the 10-digit level. 



References 



[1] ATLAS Collaboration, ATLAS Detector and Physics Performance Technical Design 
Report, Vol. 2, CERN-LHCC-99-15. 



[2] 
[3] 

[4, 



[5] 
[6] 
[7 



G. Aad et al. [The ATLAS Collaboration], larXiv:0901.05T2l [hep-ex]. 

G. L. Bayatian et al. [CMS Collaboration], J. Phys. G 34 (2007) 995. 

B. P. Kersevan and E. Richter-Was, Eur. Phys. J. C 25 (2002) 379 
|arXiv:hep-ph/0203148|; Comput. Phys. Commun. 149, 142 (2003) 
[arXiv:hep-ph/0201302 |. 

J. Cammin and M. Schumacher, ATL-PHYS-2003-024. 



V. Drollinger, T. Miiller and D. Denegri, arXiv:hep-ph/0111312 



S. Cucciarelli et al, CMS Note 2006/119; 

D. Benedetti et al, J. Phys. G 34 (2007) N221. 

J. M. Butterworth, A. R. Davison, M. Rubin and G. P. Salam, Phys. Rev. Lett. 100 
(2008) 242001 [arXiv:0802.2470 [hep-ph]]. 



[9] T. Plehn, G. P. Salam and M. Spannowsky, ' arXiv:0910.5472 [hep-ph]. 

[10] A. Bredenstein, A. Denner, S. Dittmaier and S. Pozzorini, JHEP 0808 (2008) 108 
|arXiv: 0807. 12481 [hep-ph]]. 

[11] A. Bredenstein, A. Denner, S. Dittmaier and S. Pozzorini, Phys. Rev. Lett. 103 
(2009) 012002 |arXiv:0905.0110l [hep-ph]]. 

[12] G. Bevilacqua, M. Czakon, C. G. Papadopoulos, R. Pittau and M. Worek, JHEP 
0909 (2009) 109 |arXiv:0907.4 723 [hep-ph]]. 

[13] W. Beenakker, S. Dittmaier, M. Kramer, B. Pliimper, M. Spira and P. M. Zerwas, 
Phys. Rev. Lett. 87 (2001) 201805 | arXiv:hep-ph/0107081] ; 

S. Dawson, L . H. Orr, L. Reina and D. Wackeroth, Phys. Rev. D 67 (2003) 071503 
|arXiv:hep- ph/0211438| ; 

S. Dawson, C. Jackson, L. H. Orr, L. Reina and D. Wackeroth, Phys. Rev. D 68 
(2003) 034022 [arXiv:hep-ph/0305087] . 

[14] W. Beenakker, S. Dittmaier, M. Kramer, B. Pliimper, M. Spira and P. M. Zerwas, 
Nucl. Phys. B 653 (2003) 151 [arXiv:hep-ph70211352] . 



37 



[15] S. Dittmaier, P. Uwer and S. Weinzierl, Phys. Rev. Lett. 98 (2007) 262002 
[arXiv:hep-ph/0703120] and Eur. Phys. J. C 59 (2009) 625 jarXiv:0810.0452l [hep- 
ph]]. 

[16] A. Lazopoulos, T. McElmurry, K. Melnikov and F. Petriello, Phys. Lett. B 666 (2008) 
62 |ar2£iv:Q804.2220j [hep-ph]]. 



[17] C. Buttar et al. [QCD, EW, and Higgs Working Group], arXiv:hep-ph/0604120 



[18] Z. Bern et al. [NLO Multileg Working Group]. larXiv: 0803. 04941 [hep-ph]. 

[19] A. Ferrogha, M. Passera, G. Passarino and S. Uccirati, Nucl. Phys. B 650, 162 (2003) 
(a^iv:hep-ph/0209219| . 



[20] A. Denner and S. Dittmaier, Nucl. Phys. B 658 (2003) 175 |arXiv:hep-ph/02 12259]. 



[21] A. Denner and S. Dittmaier, Nucl Phys. B 734 (2006) 62 [arXiv:hep-ph/0509141 



[22] F. del Aguila and R. Pittau, JHEP 0407, 017 (2004) |arXiv:hep-ph/0404120|. 



[23] W. T. Giele and E. W. N. Glover, JHEP 0404, 029 (2004) | arXiv:hep-ph/0402152] . 

[24] W. Giele, E. W. N. Glover and G. Zanderighi, Nucl. Phys. Proc. Suppl. 135 (2004) 
275 |arXiv: hep-ph/040701 6|; 

R. K. Ellis, W . T. Giele and G. Zanderighi, Phys. Rev. D 73 (2006) 014027 
|arXiv:hep-ph7 0508308| . 

[25] T. Binoth, J. P. G uillet, G. Heinrich, E. Pilon and C. Schubert, JHEP 0510 (2005) 
015 |arXiv:hep-ph/050426^. 

[26] T. Binoth, J. P. Guillet and G. Heinrich, JHEP 0702, 013 (2007) 
I arXiv:hep-ph/060905 4| . 

[27] T. Diakonidis, J. Fleischer, J. Gluza, K. Kajda, T. Riemann and J. B. Tausk, Phys. 
Rev. D 80 (2009) 036003 [arXiv: 0812.2134 [hep-ph]]. 

[28] T. Diakonidis, J. Fleischer, T. Riemann and J. B. Tausk, Phys. Lett. B 683 (2010) 
69 [ arXiv:0907.2115 [hep-ph]]. 

[29] Z. Bern, L. J. Dixon, D. C. Dunbar and D. A. Kosower, Nucl. Phys. B 435, 59 (1995) 
tarXiv:hep-ph/9409265j . 

[30] C. F. Berger, Z. Bern, L. J. Dixon, D. Forde and D. A. Kosower, Phys. Rev. D 74, 
036009 (2006) [arXiv:hep-ph/0604195|. 



[31] Z. Bern, L. J. Dixon and D. A. Kosower, Annals Phys. 322, 1587 (2007) 
|arXiv:0704.2798l [hep-ph]]. 

[32] C. F. Berger et al, Phys. Rev. D 78, 036003 (2008) |arXiv: 0803.4T80l [hep-ph]]. 

[33] R. Britto, F. Cachazo and B. Feng, Nucl. Phys. B 725, 275 (2005) 
[arXiv:hep-th/0412103] . 

38 



[34] A. Brandhuber, S. McNamara, B. J. Spence and G. Travaglini, JHEP 0510, Oil 
(2005) | arXiv:hep-th/0506068 |. 

[35] R. Britto, B. Feng and P. Mastrolia, Phys. Rev. D 73, 105004 (2006) 
| arXiv:hep-ph/0602178 ] . 

[36] D. Forde, Phys. Rev. D 75, 125019 (2007) |arXiv: 0704 .18351 [hep-ph]]. 

[37] S. D. Badger, E. W. N. Glover and K. Risager, JHEP 0707, 066 (2007) 



[a?Xiv:0704.3914l [hep-ph]]. 

[38] C. Ana stasiou, R. Britto, B. Feng, Z. Kunszt and P. Mastroha, JHEP 0703, 111 
(2007) [arX iv:hep-ph/0 612277:|. 

[39] W. T. Giele, Z. Kunszt and K. Melnikov, JHEP 0804, 049 (2008) [ arXiv:0801.2"237l 
[hep-ph]]. 

[40] R. K. Elhs, W. T. Giele, Z. Kunszt and K. Melnikov, Nucl. Phys. B 822 (2009) 270 
[ar Xiv:0806.3467 [hep-ph]]. 

[41] G. Ossola, C. G. Papadopoulos and R. Pittau, Nucl. Phys. B 763, 147 (2007) 
tarXiv:hep-ph/0609007| . 

[42] P. Mastroha, G. Ossola, C. G. Papadopoulos and R. Pittau, JHEP 0806, 030 (2008) 
(arXiv:0803.396l [hep-ph]]. 

[43] P. Draggiotis, M. V. Garzelli, C. G. Papadopoulos and R. Pittau, JHEP 0904 (2009) 
072 |arXiv:0903.0356l [hep-ph]]. 

[44] A. van Hameren, C. G. Papadopoulos and R. Pittau, JHEP 0909, 106 (2009) 
[arXiv:0903.4665l [hep-ph]]. 

[45] R. Keith Ellis, K. Melnikov and G. Zanderighi, Phys. Rev. D 80, 094002 (2009) 
[arXiv:0906. 14451 [hep-ph]]. 

[46] C. F. Berger et al, Phys. Rev. D 80 (2009) 074036 [ arXiv:0907. 19841 [hep-ph]]. 

[47] T. Binoth, N. Greiner, A. Guffanti, J. P. Guillet, T. Reiter and J. Renter, 
iarXiv:0910 l379l [hep-ph] . 

[48] A. Cafarella, C. G. Papadopoulos and M. Worek, Comput. Phys. Commun. 180 
(2009) 1941 [arXiv:0710.2427l [hep-ph]]. 

[49] M. Czakon, C. G. Papadopoulos and M. Worek, JHEP 0908 (2009) 085 
|arXiv:0905.0883l [hep-ph]]. 

[50] S. Catani a nd M. H. S eymour, Nucl. Phys. B 485 (1997) 291 [Erratum-ibid. B 510 
(1998) 503] |arXiv:hep-ph/9605323]. 



[51] S. Dittmaier, Nucl. Phys. B 565 (2000) 69 [arXiv:hep-ph/9904440| . 



39 



L. Phaf and S. Weinzierl, JHEP 0104 (2001) 006 |arXiv:hep-ph/0102207| . 



S. Catani, S. Dittmaier, M. H. Seymour and Z. Trocsanyi, Nucl. Phys. B 627 (2002) 
189 [arXiv:hep-ph/0201036j . 

F. A. Berends, R. P ittau and R. Kleiss, Nucl. Phys. B 424 (1994) 
308 |arXiv:hep-ph/9404313] and Comput. Phys. Commun. 85 (1995) 437 
(arXiv:hep-ph/9409326|; 

F. A. Berends, P. H. Daverveldt and R. Kleiss, Nucl. Phys. B 253 (1985) 441; 
J. Hilgart, R. Kleiss and F. Le Diberder, Comput. Phys. Commun. 75 (1993) 191. 

A. Denner, S. Dittmaier, M. Roth and D. Wackeroth, Nucl. Phys. B 560 
(1999) 33 |arXiv:hep- ph/9904472] and Comput. Phys. Commun. 153 (2003) 462 
|arXiv: hep-p h/0209330] . 



S. Dittmaier and M. Roth, Nucl. Phys. B 642 (2002) 307 | arXiv:hep-ph/0206070| . 

T. Hahn and M. Perez- Victoria, Comput. Phys. Commun. 118 (1999) 153 
|arXiv:hep-ph/9807565|; 

T. Hahn, Nucl. Phys. Proc. Suppl. 89 (2000) 231 |ar Xiv:hep-ph/0005029| . 

J. Kiiblbeck, M. Bohm and A. Denner, Comput. Phys. Commun. 60 (1990) 165; 
H. Eck and J. Kiiblbeck, Guide to FeynArts 1.0, University of Wiirzburg, 1992. 

T. Hahn, Comput. Phys. Commun. 140 (2001) 418 | arXiv:hep-ph/0012260| . 

D. B. Melrose, Nuovo Cimento XL A (1965) 181. 

G. Passarino and M. Vehman, Nucl. Phys. B 160 (1979) 151. 

G. Lei, M. Wen-Gan, H. Liang, Z. Ren- You and J. Yi, Phys. Lett. B 654, 13 (2007) 
|ar2£ivL0708.2951 [hep-ph]]. 

G. 't Hooft and M. Veltman, Nucl. Phys. B 153 (1979) 365. 

W. Beenakker and A. Denner, Nucl. Phys. B 338 (1990) 349; 

A. Denner, U. Nierste and R. Scharf, Nucl. Phys. B 367 (1991) 637. 



S. Dittmaier, Nucl. Phys. B 675 (2003) 447 | arXiv:hep-ph/0308246l . 

A. Denner, S. Dittmai er, M. Roth and L. H. Wieders, Nucl. Phys. B 724 (2005) 247 
|arXiv:hep- ph/0505042|. 

T. Stelzer and W.F. Long, Comput. Phys. Commun. 81 (1994) 357 
|arXiv:hep-p h/9401258| ; 

J. Alwali et al, JHEP 0709 (2007) 028 |ax2Civm706.2334' [hep-ph]]. 



S. Dittmaier, Phys. Rev. D 59 (1999) 016007 |arXiv:hep-ph/9805445] . 



F. A. Berends and W. T. Giele, Nucl. Phys. B 306, 759 (1988). 



40 



F. Caravaglios and M. Moretti, Phys. Lett. B 358 (1995) 332 | arXiv:hep-ph/9507237| . 

P. Draggiotis, R. H. P. Kleiss and C. G. Papadopoulos, Phys. Lett. B 439, 157 (1998) 
tarXiv:hep-ph/9807207j . 

R. Frederix, T. Gehrmann and N. Greiner, JHEP 0809 (2008) 122 [ arXiv:0808.2T28l 
[hep-ph]]. 

A. Bredenstein , S. Dittmaier and M. Roth, Eur. Phys. J. C 44 (2005) 27 
|arXiv:hep-ph/ 0506005| . 

M. Ciccohni, A. Denner and S. Dittmaier, Phys. Rev. Lett. 99 (2007) 161803 
|arXiv:0707.0381l [hep-ph]] and Phys. Rev. D 77 (2008) 013002 iarXiv:0710.4749l [hep- 
ph]]. 



T. Gleisberg et al, JHEP 0402 (2004) 056 |arXiv:hep-ph/0311263 



Tevatron Electroweak Working Group and CDF and DO Collaborations, 
a rXiv:0803.168 3 [hep-ex]. 

S. Catani, Y. L. Dokshitzer and B. R. Webber, Phys. Lett. B 285 (1992) 291. 

G. C. Blazey et al, |arXiv:hep-ex/0005012l in Proceedings of the Physics at RUN IL 
QCD and Weak Boson Physics Workshop, Batavia, Illinois, 4-6 Nov 1999, p. 47. 

[79] J. Pumplin et al, JHEP 0207 (2002) 012 |a rXiv:hep-ph/ 0201195| ; 
D. Stump et al, JHEP 0310 (2003) 046 [arXiv:hep-ph/0303013^ 



41 



